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ABSTRACT 


To  meet  the  United  States  Government  nuclear  explosion  monitoring  requirements  with  high  confidence,  the  Air 
Force  Technical  Applications  Center  needs  new  and  improved  capabilities  for  analyzing  regional  seismic, 
teleseismic,  and  infrasound  event  data.  Recently,  the  National  Nuclear  Security  Administration  has  decided  to 
investigate  three-dimensional  (3D)  modeling  in  an  effort  to  further  improve  knowledge  of  the  compressional  and 
shear-velocity  structure  as  well  as  reduce  uncertainty  and  more  accurately  detect,  locate,  and  identify  small  (body 
wave  magnitude  mb<4)  seismic  events.  For  seismically  active  areas,  inaccurate  models  can  be  corrected  using  the 
kriging  methodology;  therefore,  it  is  possible  to  detect,  locate,  and  identify  large  events  even  with  limited  resolution 
models.  This  is  not  necessarily  the  case  for  smaller  events,  however,  and  it  is  even  more  of  a  challenge  for  aseismic 
regions.  On  the  other  hand,  improving  near-regional  to  local  monitoring  demands  that  we  address  the  Earth’s 
heterogeneities  and  3D  complexities. 

Motivated  by  the  shortcomings  of  existing  single-parameter  inversion  methods  in  accurate  prediction  of  other 
geophysical  parameters,  this  research  was  mainly  focused  during  its  first  year  on  the  development  of  advanced 
multivariate  inversion  techniques  to  generate  a  realistic,  comprehensive,  and  high-resolution  3D  model  of  the 
seismic  structure  of  the  crust  and  upper  mantle  that  satisfies  multiple  independent  geophysical  datasets.  During  its 
second  year,  we  have  focused  on  the  efficient  implementation  of  the  newly  developed  technique.  Application  to 
different  areas  around  the  globe  with  different  sets  of  observations  allows  us  to  study  sensitivities,  trade-offs,  and 
possible  improvements  of  the  methodology.  We  present  3D  seismic  velocity  models  of  the  crust  and  upper  mantle 
beneath  several  regions,  resulting  from  the  simultaneous  and  joint  use  of  seismic  body-wave  arrival  times, 
surface-wave  dispersion  measurements,  and  gravity  data.  The  joint  inversion  takes  advantage  of  strengths  of 
individual  datasets  and  is  able  to  better  constrain  the  velocity  models  from  shallower  to  greater  depths.  Combining 
three  different  datasets  to  jointly  invert  for  the  velocity  structure  is  equivalent  to  a  multiple-objective  optimization 
problem.  Because  it  is  unlikely  that  the  different  “objectives”  (data  types)  would  be  optimized  by  the  same 
parameter  choices,  some  trade-off  between  the  objectives  is  needed.  The  optimum  weighting  scheme  for  different 
data  types  is  based  on  relative  uncertainties  of  individual  observations  and  their  sensitivities  to  model  parameters. 
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OBJECTIVES 

The  ultimate  goal  of  this  study  is  to  improve  our  knowledge  of  the  3D  compressional  and  shear-velocity  structure 
and  enable  us  to  reduce  uncertainty  and  more  accurately  detect,  locate,  and  identify  small  (body-wave  magnitude 
mb<4)  seismic  events,  and  therefore  improve  our  capabilities  for  nuclear  explosion  monitoring  (NEM).  This  project 
specifically  improves  seismic  monitoring  technology  through  the  development  and  application  of  advanced 
multivariate  inversion  techniques  to  generate  a  realistic,  comprehensive,  and  high-resolution  3D  model  of  the 
seismic  structure  of  the  crust  and  upper  mantle  that  satisfies  numerous  independent  geophysical  datasets. 


RESEARCH  ACCOMPLISHED 

This  is  the  second  year  of  a  three-year  project.  During  the  first  year  and  motivated  by  the  shortcomings  of  existing 
single-parameter  inversion  methods  in  accurate  prediction  of  other  geophysical  parameters  (e.g.  Julia  et  al.  2000, 
2003;  Maceira,  2006;  Maceira  and  Ammon,  2006),  we  focused  on  the  development  of  advanced  multivariate 
inversion  techniques  to  generate  a  realistic,  comprehensive,  and  high-resolution  3D  model  of  the  seismic  structure  of 
the  crust  and  upper  mantle  that  satisfies  multiple  independent  geophysical  datasets  (Maceira  et  al,  2009).  In 
particular,  we  developed  a  method  to  simultaneously  and  jointly  invert  four  different  datasets:  surface-wave 
dispersion  measurements,  teleseismic  P-wave  receiver  functions,  gravity  observations,  and  body  wave  (both  P  and  S 
phases)  travel-time  arrivals.  Surface  wave  dispersion  measurements  are  primarily  sensitive  to  seismic  shear  wave 
velocities.  Theoretically,  the  dispersion  curve  is  a  nonlinear  function  of  shear  wave  velocity  and/or  compressional 
wave  velocity  and  density  of  the  media  (Bucher  and  Smith,  1971).  However,  it  has  been  proved  (Takeuchi  and 
Saito,  1972;  Aki  and  Richards,  1980)  that  the  sensitivity  to  the  P-wave  velocity  is  significantly  smaller  than  the 
sensitivity  to  the  S-wave  velocity.  Also  the  sensitivity  function  for  the  density  is  smaller  than  the  one  for  the 
shear-wave  velocity  (Bache  et  al.,  1978;  Tanimoto,  1991).  Therefore,  shear  velocity  variations  are  usually  the  model 
parameters  in  inversion  studies  of  surface-wave  dispersion.  But,  at  shallow  depths,  it  is  difficult  to  obtain  high 
resolution  and  constrain  the  structure.  This  is  because  the  longer  the  period,  the  deeper  the  surface-wave  energy 
penetrates,  so  shorter  periods  are  primarily  sensitive  to  upper  crustal  structures.  But  short  periods  are  difficult  to 
measure,  especially  in  tectonically  and  geologically  complex  areas.  On  the  other  hand,  gravity  inversions  have  the 
greatest  resolving  power  at  shallow  depths  since  gravity  anomalies  decrease  in  amplitude  and  increase  in  wavelength 
with  increasing  depth.  And  gravity  measurements  supply  constraints  on  rock  density  variations.  In  addition, 
surface-wave  dispersion  measurements  are  primarily  sensitive  to  vertical  shear-wave  velocity  averages;  while  body 
wave  receiver  functions  are  sensitive  to  shear-wave  velocity  contrasts  and  vertical  travel-times.  Addition  of  the 
seismic  travel-time  data  helps  to  constrain  the  shear  wave  velocities  both  vertically  and  horizontally  in  the  model 
cells  crossed  by  the  ray  paths.  Thus  by  combining  these  four  complimentary  datasets  into  a  single  inversion,  we  can 
obtain  a  self-consistent  3D  seismic  velocity-density  model  with  increased  resolution  of  shallow  geologic  structures. 
Our  final  algorithm  is  a  modification  of  the  Maceira  and  Ammon  (2009)  joint-inversion  code,  in  combination  with 
the  regional  version  of  the  double-difference  (DD)  tomography  program  tomoDD  (Zhang  and  Thurber,  2003,  2006), 
with  a  fast  LSQR  solver  operating  on  the  gridded  values  jointly. 

During  this  second  year,  we  have  focused  on  the  efficient  implementation  of  this  newly  developed  technique. 
Application  to  different  areas  around  the  globe  with  different  sets  of  observations  allows  us  to  study  sensitivities, 
trade-offs,  and  possible  improvements  of  the  methodology.  We  present  here  3D  seismic  velocity  models  of  the  crust 
and  upper  mantle  beneath  several  regions,  resulting  from  the  simultaneous  and  joint  use  of  seismic  body-wave 
arrival  times,  surface-wave  dispersion  measurements,  and  gravity  data. 


East  Africa  Rift  System 

Knowledge  of  crustal  and  upper  mantle  structure  is  of  importance  for  understanding  East  Africa’s  geodynamic 
evolution  and  for  addressing  broader  questions  about  the  causes  of  continental  breakup.  Though  recent 
investigations  have  yielded  improved  characterizations  of  the  rift  zone,  questions  remain  concerning  the  distribution 
of  velocity  anomalies  (e.g.,  Julia  et  al.,  2005),  extent  of  lithospheric  thinning  (e.g.,  Weeraratne  et  al.,  2003),  and 
temporal  variation  in  rift  evolution  (e.g.,  Ebinger,  1989).  Key  to  resolving  these  questions  are  better  constrained 
seismic  models.  We  decided  to  test  our  joint  inversion  methodology  in  this  region  due  to  its  scientific  significance 
and  because  of  a  more  practical  reason:  the  unprecedented  opportunity  of  getting  a  handle  on  on-land  Bouguer 
gravity  observations.  Benefits  of  our  joint  inversion  approach  appear  pronounced  when  working  with  regions  of 
strong  lateral  contrast  as  found  in  central  Asia  (Maceira  and  Ammon,  2009).  In  applying  the  joint-inversion 
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technique  to  East  Africa,  we  solve  for  velocity  structure  in  an  area  with  less  lateral  heterogeneity  but  great  tectonic 
complexity.  To  increase  the  effectiveness  of  the  technique  in  this  region,  we  explore  gravity-filtering  methods  and 
test  different  velocity-density  relations  (the  relationship  between  the  independent  datasets  is  always  one  of  the 

difficulties  with  simultaneous  joint  inversions). 

The  area  for  the  inversion  spans  the  broad  uplifted  region  from 
Ethiopia  at  one  end  to  Kenya,  Uganda,  and  Tanzania  at  the 
other  (Figure  1).  Near  the  northern  boundary  of  our  study  area, 
the  Main  Ethiopia  Rift  meets  the  incipient  Red  Sea  and  Gulf  of 
Aden  spreading  ridges.  At  the  opposite  end,  the  rift  system 
splits  into  distinct  western  and  eastern  branches,  which  largely 
sidestep  the  Archean  Tanzania  craton.  Recent  inversions  of 
East  Africa  have  employed  body  waves  (e.g.,  Benoit  et  al., 
2006),  surface  waves  (e.g.,  Weeraratne  et  al.,  2003),  receiver 
functions,  or  some  combination  of  these  (e.g.,  Julia  et  al., 

2005;  Keranen  et  al.,  2009).  Although  useful  comparisons  can 
be  drawn  between  the  Ethiopian  and  Tanzanian  portions  of  the 
rift  system,  most  tomographic  studies  to  date  have  focused 
exclusively  on  one  section  or  the  other.  The  current  inversion, 
in  contrast,  is  carried  out  over  a  wider  area  than  most  previous 
studies,  allowing  straightforward  comparison  between  these 
two  distinct  portions  of  the  rift  system. 

Fundamental-mode  Rayleigh  wave  group  velocity  estimates 
with  periods  from  7  s  to  150  s  were  obtained  from  Pasyanos 
and  Nyblade  (2007)  for  the  inversion.  Though  less  detailed 
than  images  from  local  seismic  arrays  (e.g.,  Prodehl  et  al., 
1997),  these  estimates  span  a  broader  spatial  and  period  range. 
Gravity  data  for  the  inversion  were  derived  from  the  Gravity 
Recovery  and  Climate  Experiment  (GRACE)  (Tapley  et  al., 
2005).  Comparison  of  GRACE-converted  Bouguer  anomalies 
with  values  from  a  proprietary  gravity  compilation  showed 
close  agreement  between  the  two  datasets;  in  light  of  this 
similarity,  we  chose  the  GRACE  data  because  of  their  more 
uniform  coverage. 

We  developed  and  implemented  a  method  to  increase  the 
usefulness  of  gravity  data  by  filtering  the  Bouguer  anomaly 
map.  Though  commonly  applied  in  gravity  forward  modeling 
(e.g.,  Simiyu  and  Keller,  1997)  such  techniques  have  not  to  our 
knowledge  been  used  in  previous  joint  inversion  studies  (e.g.,  Lees  and  VanDecar,  1991;  Zeyen  and  Achauer  1997; 
Tiber!  et  al.,  2003;  Maceira  and  Ammon,  2009).  Figure  2a  suggests  possible  benefits  of  filtering.  In  this  figure, 
results  from  the  inversion  using  group  velocities  only  (solid  line)  provide  the  starting  model  for  the  joint  inversion 
(dashed  line).  Comparison  of  the  two  profiles  suggests  that  the  addition  of  unfiltered  gravity  data  contributes  little  in 
the  way  of  distinguishing  between  features  at  different  depths.  Rather  than  improving  resolution  at  shallow  depths  as 
desired,  features  in  the  unfiltered  gravity  data  are  smeared  into  the  mantle;  moreover,  within  the  mantle  the  two 
profiles  are  offset  by  a  roughly  constant  amount,  suggesting  that  gravity  data  provide  inadequate  constraints  on 
mantle  anomalies.  Examination  of  additional  profiles  shows  that  these  problems,  likely  exacerbated  by  East  Africa’s 
pronounced  long-wavelength  Bouguer  components  (Girdler,  1975),  are  pervasive  throughout  the  study  area.  In 
response,  we  employ  filtering  methods  to  remove  the  long-wavelength  gravity  components,  reducing  the  tendency 
of  crustal  gravity  features  to  be  smeared  into  the  mantle.  To  balance  the  change  to  the  observed  gravity  values 
brought  on  by  filtering,  we  reduce  the  contribution  to  the  predicted  gravity  values  from  blocks  below  60  km, 
corresponding  to  the  portion  of  the  Bouguer  anomaly  signal  we  attempt  to  suppress.  Finally,  while  allowing  velocity 
structure  above  60  km  to  vary  freely,  below  60  km  we  induce  agreement  between  the  joint  inversion  results  and  the 
initial  model  by  imposing  stronger  a  priori  model  constraints.  In  the  application  of  the  described  technique, 
problems  of  depth  correlation  associated  with  filtering  make  our  choice  of  depth  cutoff  (60  km)  somewhat 
subjective.  This  disadvantage,  however,  is  likely  outweighed  by  the  reduction  in  smearing  of  the  gravity  signal 
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Figure  1.  We  distinguish  between  Ethiopia  and 
Tanzania  as  distinct  portions  of  the 
larger  rift  system.  The  Ethiopian 
volcanic  province  (after  Keranen  et 
al.,  2009)  spans  the  northern  portion 
of  the  rift  system.  The  Tanzania 
craton  (after  Nyblade  and  Brazier, 
2002)  marks  the  southern  portion  of 
the  rift  system  (rift  lies  for  the  most 
part  along  the  edges  of  the  craton). 
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suggested  by  Figure  2b.  Although  filtering  removes  potentially  useful  information  on  mantle  structure,  the 
remaining  short- wavelength  signal  can  be  assigned  with  greater  reliability  within  the  crust,  avoiding  the  mutually 
degrading  effects  of  smearing  between  crust  and  mantle.  To  remove  the  long-wavelength  components  from  the 
Bouguer  gravity  map,  we  follow  Tessema  and  Antoine  (2004),  who  use  an  upward  continuation  method  and 
demonstrate  correlation  with  crustal  geology. 


40  km 


60  km 


80  km 


Figure  2.  One-dimensional  profiles  from  a  representative  cell 
in  our  model.  Comparison  shows  results  from  an 
inversion  of  group  velocities  only  (solid  lines)  and 
from  joint  inversion  (dashed  lines),  (a)  Without 
gravity  filtering,  results  point  to  an  inability  to 
distinguish  between  density  anomalies  at  depth. 

(b)  With  gravity  filtering,  the  two  profiles  agree  less 
closely  at  crustal  depths  and  more  closely  over 
mantle  depths. 

Previous  inversions  have  incorporated  gravity  data  by  formulating 
the  forward  problem  in  two  steps:  (1)  velocity  variations  are  mapped 
to  density  variations  and  (2)  density  variations  are  mapped  to  gravity 
anomalies.  While  the  second  mapping  is  mathematically 
straightforward,  the  first  necessitates  use  of  empirical  formulae. 
Perhaps  because  empirical  mappings  between  density  and  P-wave 
velocity  (e.g.,  Birch’s  law)  are  more  familiar  than  mappings  between  density  and  S-wave  velocity,  previous  studies 
have  relied  on  the  former  together  with  ad  hoc  scaling  to  obtain  S-wave  velocities.  Here,  we  eliminate  the  need  for 
such  scaling  by  adopting  a  direct  mapping  between  density  and  S-wave  velocity.  Given  the  variability  in  Poisson’s 
ratio  for  rift  zones  (e.g.,  Dugda  et  al.,  2005),  this  new  mapping  provides  the  proper  framework  for  allowing  spatial 
Vp/Vs  variability,  which  we  envision  developing  in  future  inversions. 

Figure  3  shows  the  3D  S-wave  velocity 
model  obtained  from  the  joint  inversion. 
Two  sets  of  velocity  anomalies,  one 
beneath  Ethiopia  and  the  other  beneath  the 
Tanzania  craton,  stand  out  among  the  most 
prominent  features.  Beneath  Ethiopia  a 
distinct  low-velocity  anomaly  is  observed 
over  a  depth  of  30-150  km.  Though  the 
existence  of  lower-than-average 
shear-wave  velocities  in  the  uppermost 
mantle  beneath  Ethiopia  matches  a  variety 
of  previous  results  (e.g.,  Keranen  et  al., 
2009),  controversy  exists  regarding 
velocity  structure  at  greater  depths.  Such 
disagreement  may  in  part  reflect 
differences  in  methodology.  While  surface 
wave  studies  (e.g.,  Ritsema  and  van  Heijst, 
2000)  have  often  failed  to  resolve  low 
velocities  beneath  Ethiopia  at  depths 
greater  than  170  km,  body- wave 
investigations  (e.g.,  Benoit  et  al.,  2006) 
have  imaged  low  velocities  extending  as 
deep  as  400  km.  As  a  result  of  upward 
trajectories  associated  with  incoming  rays, 
body- wave  studies  may  tend  to 
underestimate  velocity  deviations  and 
smear  velocity  structure  (see  Weeraratne  et 

al.,  2003,  for  analysis  based  on  observations  from  East  Africa).  From  these  considerations,  it  seems  likely  that  the 
true  vertical  extent  of  the  low- velocity  anomaly  beneath  Ethiopia  lies  between  the  shallower  and  deeper  values 


Figure  3.  Horizontal  cross  sections  at  various  depths  through  the 
3D  shear  velocity  model  (velocity  units:  km/s). 
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predicted  using  surface-wave  and  body- wave  methodologies,  respectively.  Besides  the  low- velocity  anomaly 
beneath  Ethiopia,  prominent  velocity  excursions  also  occur  beneath  Tanzania.  As  shown  in  Figure  3,  we  resolve 
high  velocities  below  the  craton  over  a  depth  of  50-130  km.  Contiguous  to  this 

high-velocity  zone,  a  low-velocity  anomaly  emerges  at  a  depth  of  140  km.  Although  the  upper  part  of  this 
low-velocity  anomaly  is  well  constrained  by  the  long-period  Rayleigh  wave  group  velocities  used  in  the  inversion, 
sensitivity  kernels  suggest  resolution  drops  off  below  175  km.  This  juxtaposition  of  high  and  low  shear-wave  speeds 
appears  consistent  with  the  hypothesis,  discussed  in  detail  by  Weeraratne  et  al.  (2003),  of  a  hot  mantle  plume 
impeded  by  cool,  stable,  overlying  cratonic  lithosphere.  In  the  following  particulars,  we  obtain  close  agreement  with 
Weeraratne  et  al.  (2003):  (1)  the  low-velocity  anomaly  beneath  the  craton  originates  at  around  140  km,  (2)  the 
anomaly  lies  squarely  beneath  the  craton,  and  (3)  deviations  from  background  shear  velocities  reach  as  much  as  5% 
over  the  depths  displayed  in  Figure  3.  Additionally,  our  results  allow  comparison  between  the  rift  structures  of 
Ethiopia  and  Tanzania.  In  obtaining  data  from  Pasyanos  and  Nyblade  (2007),  we  use  group  velocities  derived  not 
only  from  local  stations  and  events,  but  also  from  stations  and  events  distributed  across  surrounding  tectonic  plates. 
Though  the  resulting  continental-scale  maps  possess  less  detail  than  local-scale  group-velocity  maps,  their  wider 
spatial  coverage  allows  straightforward  comparison  between  distinct  portions  of  the  rift  system.  As  a  result,  we  find 
that  uppermost  mantle  shear  velocities  beneath  Ethiopia  appear  much  slower  than  those  beneath  Tanzania.  While  the 
presence  of  shallow  low  velocities  beneath  Ethiopia  can  be  explained  in  terms  of  lithospheric  erosion  or 
magma-assisted  rifting  (e.g.,  Dugda  et  al.,  2007;  Ebinger  and  Casey,  2001),  the  absence  of  shallow  low  velocities 
beneath  Tanzania  suggests  an  alternate  mechanism  of  extension.  Such  differences  may  stem  in  part  from  inherited 
rheology,  as  related  for  example  to  the  high-velocity  zone  observed  beneath  Tanzania.  Finally,  though  a  common 
origin  below  175  km  cannot  be  ruled  out,  no  evidence  is  found  to  suggest  the  various  low- velocity  anomalies  in 
Figure  3  merge  continuously  at  depth.  Indeed,  low  velocities  beneath  Ethiopia  appear  diffuse  and  diminished  by 
1 60  km,  suggesting  the  existence  of  a  seismic  velocity  discontinuity  separating  the  low- velocity  anomaly  beneath 
Ethiopia  from  the  one  below  Tanzania. 
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(a)  Observations 
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Figure  4.  Fit  to  data  from  inversion  of  surface-wave  dispersion  only  and  from  joint  inversion,  (a)  Top:  Group 
velocities  from  a  representative  cell  in  our  model.  Bottom:  Filtered  Bouguer  anomalies. 

(b)  Top:  Group  velocity  fit  obtained  from  the  inversion  of  group  velocities  only.  Bottom:  Gravity  fit 
from  the  same  inversion,  (c)  Top:  Group  velocity  fit  obtained  from  joint  inversion.  Bottom:  Gravity 
fit  from  joint  inversion. 

In  the  inversion  carried  out  in  central  Asia  by  Maceira  and  Ammon  (2009),  addition  of  gravity  data  dramatically 
improved  the  fit  to  the  Bouguer  anomalies  without  significantly  degrading  the  fit  to  the  group  velocities.  Figure  4 
demonstrates  this  result  for  the  current  study  area  as  well.  Although  it  is  well  known  that  problems  of  nonuniqueness 
make  gravity  data  easier  to  match  than  seismic  data,  several  observations  provide  confidence  in  our  methodology’s 
robustness.  These  include  the  simultaneous  fit  to  both  datasets  shown  in  Figure  4  as  well  as  other  changes  resulting 
from  the  addition  of  gravity  data.  Compared  with  results  from  the  inversion  of  group  velocities  only,  the 
joint-inversion  methodology  provides  increased  effectiveness  in  capturing  Moho  depth  and  improved  resolution  of 
uppermost  mantle  detail,  including  sharper  delineation  of  the  Tanzania  craton  in  the  vertical  cross  sections  and  in  the 
horizontal  cross  section  at  60  km.  The  resolved  extent  of  the  high-velocity  cratonic  region  accords  well  with 
previous  tomographic  images  (e.g..  Figure  1 1  of  Weereratne  et  al.,  2003)  and  with  geodynamic  models  suggesting 
strain  localization  in  zones  of  weakness  surrounding  the  craton  (e.g.,  Nyblade  and  Brazier,  2002). 


Utah  Geothermal  Field 

The  Cove  Fort-Sulphurdale  geothermal  area  is  located  in  the  transition  zone  between  the  Basin  and  Range  to  the 
west  and  the  Colorado  Plateau  to  the  east  (Figure  5).  We  have  collected  various  geophysical  data  around  the 
geothermal  field,  including  heat  flow  (Henrikson  and  Chapman,  2002),  gravity  (Pan  American  Center  for  Earth  and 
Environmental  Studies  (PACES)  available  at  http://gis.utep.edu),  seismic  surface-wave  phase  and  group  velocity 


102 


2010  Monitoring  Research  Review:  Ground-Based  Nuclear  Explosion  Monitoring  Technologies 


maps  (Yang  et  al,  2008),  and  seismic  body-wave  arrival  times  that  were  assembled  from  seismic  waveforms 
recorded  by  the  University  of  Utah  Seismograph  Stations  (UUSS)  regional  network  for  the  past  7  years  and  the 
recent  EArthscope/US  Array  phase  data.  All  of  these  geophysical  datasets  have  different  strengths  for  characterizing 
subsurface  structures  and  properties.  Combining  these  data  through  a  coordinated  analysis  and,  when  possible,  by 
joint  inversion  provides  a  detailed  model  of  the  Cove  Fort  geothermal  region. 


Figure  5.  Utah  simpliBed  geological  map  (after  Stokes,  1977;  Hecker,  1993) 


Various  geophysical  datasets  indicate  that  beneath  the  Cove  Fort-Sulphurdale  geothermal  resource  there  is  a  strong 
anomaly  of  low  seismic  velocity,  low  gravity,  and  high  electrical  conductivity  that  correlates  with  the  high  surface 
heat  flow.  This  suggests  that  there  is  a  heat  source  in  the  crust  beneath  the  geothermal  field.  We  collected  first 
P-arrival  data  from  more  than  6500  earthquakes  in  the  Utah  region.  Each  event  has  at  least  six  arrivals  for  reliably 
determining  its  location.  We  apply  the  DD  seismic  tomography  method  (Zhang  and  Thurber,  2003)  to 
simultaneously  determine  an  initial  velocity  structure  and  earthquake  locations.  The  finite-difference  travel-time 
calculation  method  of  Podvin  and  Eecomte  (1991)  is  used  to  calculate  travel  times  between  events  and  stations.  On 
the  preliminary  regional  seismic  velocity  map  computed  this  way,  we  can  also  identify  some  other  low- velocity 
anomalies,  indicating  other  potential  geothermal  prospects.  We  decided  to  apply  our  joint  inversion  methodology  to 
produce  a  better-constrained  velocity  structure  of  the  Utah  area,  which  will  be  very  helpful  for  characterizing  and 
exploring  existing  and  potential  geothermal  reservoirs  in  the  area. 
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Figure  6.  Shear-wave  velocity  model  at  constant  depth  slices.  The  depth  of  each  image  is  shown  at  the 
top  of  each  map  (velocity  units:  km/s). 
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Relative  w&ighting  between 
gravity  and  surface  wava 


Relative  weighting  between 
travel  time  and  gravity 


Figure  6  shows  the  results  from  using  the  surface  and 
gravity  data  to  jointly  determine  the  shear- velocity 
(Vs)  model.  Combining  both  datasets  to  jointly  invert 
for  the  velocity  structure  is  equivalent  to  a 
multiple-objective  optimization  problem.  Because  it  is 
unlikely  that  the  different  “objectives”  (data  types) 
would  be  optimized  by  the  same  parameter  choices, 
some  trade-off  between  the  objectives  is  needed.  The 
optimum  weighting  scheme  for  different  data  types  is 
based  on  relative  uncertainties  of  individual 
observations  and  their  sensitivities  to  model 
parameters.  Here  we  determine  their  relative 
weightings  through  a  trade-off  analysis  of 
surface-wave  residual  and  gravity  data  residual 
(Figure  7).  As  a  result,  the  model  optimally  fits  both 
the  surface  wave  and  gravity  data.  In  this  study,  the 
surface  wave  is  weighted  50  times  of  the  gravity  data 
(1:0.02).  Compared  to  the  inversion  using  only  the  surface  wave  data,  at  shallow  depths  (e.g.,  1  and  5  km),  the  Vs 
model  better  delineates  the  low-velocity  anomaly.  The  first-order  horizontal  and  vertical  smooth  weighting  are  also 
applied  to  constrain  the  model  during  the  joint  inversion.  After  three  iterations  of  joint  inversion,  the  surface-wave 
data  residual  is  close  to  zero  and  the  gravity  data  residual  decreases  -90%. 

Preliminary  joint  inversion  results  using  seismic  travel  time  and  gravity  data  (Figure  8)  indicate  strong  low  velocity 
anomalies  in  middle  crust  beneath  some  known  geothermal  sites  in  Utah. 


Figure  7.  Trade-off  curves  between  surface-wave  and 
gravity  data  residuals  (left)  and  between 
travel  time  and  gravity  data  residuals 
(right). 
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Figure  8.  Compressional  wave  velocity  model  at  constant  depth  slices.  The  depth  of  each  image  is  shown  at 
the  top  of  each  map  (velocity  units:  km/s). 
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CONCLUSIONS  AND  RECOMMENDATIONS 


This  is  the  second  year  of  a  three-year  project  to  map  the  3D  seismic  structure  of  the  crust  and  upper  mantle  using 
simultaneous  joint  inversion  of  surface-wave  dispersion,  gravity,  receiver  function,  and  travel-time  observations. 
Geophysical  models  play  an  important  role  in  ground-based  nuclear  explosion  monitoring  (GNEM).  To  more 
confidently  and  accurately  detect,  locate,  and  identify  small  seismic  events,  better  high-resolution  3D  structural 
models  are  needed.  Therefore,  the  ongoing  research  directly  addresses  this  challenge,  and  our  results  are  also  of 
fundamental  importance  for  understanding  the  geodynamic  evolution  and  formation  of  continents,  as  well  as  the 
processes  acting  within  and  on  the  continental  lithosphere. 

We  are  now  focusing  on  developing  3D  models  for  different  areas  for  which  we  have  different  sets  of  observations. 
This  approach  allows  us  to  study  sensitivities,  trade-offs,  and  possible  improvements  of  our  methodology.  We  have 
presented  here  models  obtained  for  the  East  Africa  Rift  System  and  the  Utah  geothermal  area  (more  areas  of  interest 
to  the  GNEM  program  will  be  shown  at  the  MRR  2010  meeting).  By  tailoring  our  methodology  to  these  two 
different  areas,  we  have  learned  that  besides  enhancing  resolution  at  short  wavelengths,  use  of  filtered  gravity 
anomalies  may  help  distinguish  between  anomalies  at  different  depths.  We  have  also  tested  different  relations 
between  seismic  velocities  and  density  and  concluded  that  a  relation  where  no  a  priori  Poisson’s  ratio  assumptions 
are  needed  (e.g.,  D.  G.  Harkrider’s)  is  more  beneficial.  Future  work  will  involve  validation  of  these  3D  models 
through  the  use  of  the  spectral  element  method  and  full  waveform  comparisons. 
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